Effective potentials for heavy quarks above deconfinement 

K. Duslin^ and C. Youn^ 

Department of Physics & Astronomy, State University of New York, 
Stony Brook, NY 11794-3800, U.S.A. 
(Dated: October 14, 2009) 

Abstract 

Using methods previously developed by Kelbg and others for creating effective potentials for 
electron-ion plasmas, we investigate quarkonium potentials above deconfinement. Using results 
for the internal energy of a static quark-antiquark (QQ) pair, a set of effective potentials are 
constructed that take into account quantum effects and symmetry conditions. Bound state effects 
are explicitly included in order to account for the strongly coupled nature of the plasma. It is hoped 
that these effective potentials will be useful in simulations of heavy quarks or plasma quasiparticles 
when the dynamics is treated classically. 
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I. INTRODUCTION 



It has recently been noted that the matter produced in heavy ion colhsions at RHIC 
cannot be weakly coupled but instead behaves as a good liquid [H [2] . The evidence for the 
sQGP is large but mainly consists of the following two points: 1. the observed collective 
flows at RHIC can be explained by hydrodynamics showing that the dissipative lengths are 
very short and 2. binary bound states are seen to exist in lattice simulations above Tc and 
are also predicted 13] using lattice interparticle potentials. 

Since a complete quantum description tends to rely on perturbative methods other ap- 
proaches have to be adopted in order to perform calculations at strong coupling. One 
approach as was discussed in [H El E] , is to model the strongly interacting quark and gluon 
quasiparticles as a classical non-relativistic colored Coulomb gas. This model is analyzed 
using molecular dynamics (MD) simulations in which real time correlators can be extracted. 
The interaction between quarks was taken as a Coulomb potential with a strong repulsive 
core in order to mock up quantum-mechanical effects. 

A second example is a model treating charmonium in the sQGP [7|, where initial ensem- 
bles of quark pairs from pQCD events are evolved classically according to a Langevin model 
with an additional interaction potential included between quark-antiquark pairs. The in- 
teraction potential was taken from lattice cc internal energy data with quantum-mechanical 
effects mimicked by simply turning off the potential below the approximate Bohr radius. 

If these models are to be refined it is necessary that we consider quantum-mechanical 
effects more carefully. The construction of an effective interaction which mimics the effects 
of quantum dynamics was first performed by Kelbg [10] . The goal of this work is to generate 
classical effective qq potentials which take into account quantum-mechanical effects, first the 
dynamics and then the symmetry or anti-symmetry which must exist in the wavefunction 
for the qq pair. We want to find an effective potential Veff{r, T) whose classical Boltzmann 
factor yields the diagonal term of the quantum — mechanical density matrix: 



For example, if this corrected interaction were used in the modeling of charmonium discussed 
in the previous paragraph, the final equilibrium distribution obtained by the Fokker-Planck 
evolution would be the correct quantum-mechanical distribution. This would obviously 




(1) 
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improve estimates for the intermediate states as well. These effective interactions would 
also have use in the MD simulations discussed earlier. 

In this work we calculate effective quark-antiquark potentials from the internal energy of 
static QQ pairs calculated in [8] based on lattice calculations of the QQ free energy in [9]. 
Section [njof this paper summarizes a number of previous works [HI [12], [131 [H], including the 
definition of the Slater sum and its use in creating effective potentials. The original results 
of this paper start in section III, where the methods presented in section [II] are applied to 



potentials extracted from lattice data. 



II. CREATING EFFECTIVE POTENTIALS 



A. The Slater Sum 



Anticipating applications of this work to many-particle systems, we define the A^-particle 
Slater sum: 

^W(fi,...,r-W) = iV!A3^^|vI/,(rl,...,r1v)re-^^« (2) 



where are the N-particle energy eigenf unctions and A = y 2i!mkT thermal wave- 

length. Integrating the Slater sum in coordinate space yields the partition function: 

Z = I (c?n...rfr-W)^(^)(n,...,r-V) (3) 

In the case of ideal particles using Boltzmann statistics the partition function is given as 
Zmb = N\\3N ■ The corresponding free energy can then be separated in the following way: 

F = -kBTlnZ = Fmb - ksT In (4) 

^MB 

Substituting eq. |3] for Z in the above equation yields: 

F = Fmb - ksTln ^ J {dn...dfN)S^''\n, ...,f^) (5) 
This expression should now be compared to the classical counterpart of the free energy: 

Fd = Fmb - ^sTln ^ J {dn, drW)e-^S«<. (6) 
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In the classical limit one obtains: 

4^)(n,...,r-W)=e-^^-.^- (7) 

and near the classical limit one can include multi-body interactions due to the quantum 
correlations between particles: 

S^^\fi, ...,rW) = e-^^»<^-"'^-'^^»<:'<'="^^-^-+- (8) 

where Uij and Uijk are effective two- and three-body interactions respectively and can be 
expressed in terms of two- and three-particle Slater sums as: 



u. 



Keeping only the first term, Uij, in the series for S*^^-* accounts for quantum effects up to 
first order in the coupling parameter, T = where Ze,aws,T are respectively the ion 

charge, the Wigner-Seitz radius aws = (3/47™)^/'^ and the temperature. The second term 
Uijk introduces corrections to the effective potential of order F^. 



B. Effective potentials in the binary approximation 

The previous section motivates us to take the binary approximation: 

5(A^) ~e-^^«<^ "'^- (10) 

where Uij is defined as in eq. [9] This approximation is exact for a 2-body systems (for 
example the model discussed in section [T] where a heavy quark only interacts with its diagonal 
partner). For the remainder of this paper, we now focus on an explicit calculation of Uij. 
Our method involves examining the two-particle density matrix: 

p(ri, f2, f/, f^) = J2 (^1, r2)^i(f/, r2')e-^^» (11) 

i 

where the diagonal term is proportional to the needed Slater sum. The two- 
particle density matrix can then be factored into center-of-mass and relative components, 
p(ri,f2,f/,f2) = p{R,R')p{f,f'). 
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The two-body effective potential can then be found in the hmit as R' ^ R and f' f. 
Since p{R,R) = 1 one is left with solving for p{r,r). 

Acting on p(r', r') with the one-particle Hamiltonian^, H = — 2^V^ + V{r) where = 
mj~^ + ^ is the reduced mass of the two particle system and r is the relative coordinate 
it can be shown that p satisfies: 

H{r)p{f,f^,(3) = -^ (12) 

The boundary condition on the density matrix is given from the fact that {"^i} is a 
complete basis: p(r, r /? = 0) = (5(r — r ') 



We use the same method of solution on eq. 12 which is given in [TT]. The solution to the 



free-particle density matrix (V=0 in eq. 12) is given as 



1 — |r — f 



/|2 



where ^ = y 2^ reduced mass of the two-particle system. 

We then define the effective potential as in the introduction: 

p(f, f', /3) = poe"'^"^"''"''^^ (14) 

The effective two body potential Uij of interest is simply the diagonal part of u{r, t\ (3) when 

Substituting the above form for the density matrix 14 into the Bloch equation 12 one 
finds that u{r,r',P) satisfies: 

— /3V\ - —p\Vu? - (f - f') ■ Vm + V(r) = muld(3) + u (15) 
2m 2m 



In the limit that m, VrW, or /? is small, eq. 15 may be linearized: 

^/3V\i - (f - r') ■ Vmi + ^(f) = (3{dui/d(3) + mi (16) 

the solution of which we call Ui{f,r\f3). As has been emphasized in [11], this condition is 
far more permissive than the typical condition from perturbation theory that u be small. 
This equation can be solved exactly: 

«i(F,f')= [ dfiG{r,f',n,P)V{fi) (17) 



^ From here on we use standard high-energy units where h — c — 1 
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where the Green function G is given by: 

ti\2 



G(f,r',ri,/3) 



m\r — r 



exp( 



2/3 













r — r*! 







exp( 



m(|r — r*i| + \ f' — fi\ 
2^ 



(18) 

For a spherically symmetric potential V{r), the angular integration in the above equation 
can be performed and the final result for the diagonal term is given as ui{f = f' ,(3) is: 



ui{f,j3) 



X 



ri 
r 



V{ri) [erf( 



A 



erf(- 



ril 



A 



]dri 



(19) 



Equation [19] will be the starting point for our future analysis of lattice potentials. It can 
be shown that the same result can be found by solving for the Slater sum S"*^^^ explicitly 
using plane waves for the wavefunction and solving perturbatively in the couphng constant. 



1. Bound States 



It has already been shown in ^T] that the linearized approximation which yielded eq. 16 is 
satisfactory assuming the contribution from bound states is small. The term dropped from 
the linearized equation, |/3^(Vm)^, becomes significant when the gradient of the effective 
potential is on the order of the temperature. This is exactly the condition leading to large 
contributions from the lowest bound states to the density matrix. In order to take into 
account the effect of bound states explicitly we follow [12j, where the bound and free states 
are considered as two separate contributions to the Slater sum. Both the partition function 
and the two-particle Slater sum can be separated as 

5(2)(r,/5) = {1-P')S^^^ + P'S^^\ (20) 

where the operator P' projects out or removes the free component. The terms on the 
right-hand side correspond to the free and bound state contribution to the Slater sum, 
respectively. 

There are various ways of defining the projection operator P', so that it separates the 

(2) 

Slater sum into a part dominated by the free contribution (call it Sj ) and a part dominated 
by the bound state contribution. In principle, these are merely different conventions, however 

(2) 

since we will identify Sjr with the Slater sum we already calculated, some conventions are 
better than others. Explicitly, one can use the Riewe-Rompe convention which acts as a 
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sharp cutoff between the bound and continuum states. More specifically, the projection 
operators acting on the two particle Slater sum yield 

S^'\r,P) = Sf\r,P) + J2 \^^ir)\'e-f'''\ (21) 

Ei<e* 

where e* ~ T. Using this convention, we would be adding a term to the Slater sum even 
where (3 is small, where our calculation needs no correction. In order to circumvent this, 
we use the Brillouin-Planck-Larkin (BPL) convention since it has a continuous transition at 
the continuum edge, 

S^'\r,(3) = Sf\r,P) + \'^^{r)\'[e-^'^^ - I + PE,]. (22) 

i 

The advantage of this convention can be seen by examining small and large (3. When /3 
is small, the coefficient e~^^^ — 1 + jSEi ^ (1 — l3Ei) — 1 + jSEi ^ 0, and no significant 
contribution is added to the free Slater sum. When (3 is large, e~^^^ — 1 + (3Ei ^ e~^^' {Ei 
is negative), making the bound state part the dominant contribution to the Slater sum, as 
it should be at low temperatures. 

In the above expressions S^j^\r,(3) is the free-particle Slater sum and can be found by 



substituting the result from eq. 19 into eq. p\ The second term includes a sum over bound 



states with wavefunctions given as solutions to the Schrodinger equation. 
2. Symmetry Considerations 

We make one final consideration in creating effective potentials: symmetry conditions for 
identical particles. In the case of the identical fermions, the two-particle density matrix, 
written as an imaginary-time propagator, becomes 

pin , fa, f/, f^,(3) = (fi , fa I - (fa, fi I ) exp{-f3H) ( |f/, f^) -\f^, f/) ) (23) 

For our situation where two particles interact according to a potential only dependent on 
|f|, their relative separation, this formula dramatically simplifies: 

p(fi,f2,f/,f2,/3) = pabs{R,R\(3)[prei{r,r\l3) - Prei{r, -r\ fi)] (24) 

Again, we will be interested in the limit f ' = f. We may calculate explicitly the off-diagonal 



term of the density matrix needed in eq. 24, or we may follow the approach in [13] and 



approximate the off-diagonal term using only the off diagonal free particle density matrix 

Prei{r, -f, (3) ^ exp(-2m|r1 V/5)Pre«(^ ^> /5) (25) 
The final form for the correctly symmetrized density matrix psym is given as: 

Ps,™(r,/3) = ^[l±cxp(-2m|rlV/3)]e-^"^(^^ (26) 
where a minus sign is used for fermions and a plus sign is for bosons. 

C. Coulomb Potential 

We now apply the methods outlined in the previous section to the case of the Coulomb 



potential. We first note that equation 19 can be integrated exactly in the case of a Coulomb 
potential, V{r) = — 1/r with the result: 

«i(r,/3) = ^(1 - e-^VA^ + 4^ [1 - erf(^)]) (27) 

In figure [T] we show the results of including bound states for the case of the Coulomb 
potential. The dotted green and blue curves show the results for ui{r) for values of /3 = 10 
and 1 respectively. The corresponding solid curves show the results of including the first 
lowest bound state in the sum in eqn |22j The solid red curve shows the result dX (3 = oo 
where the result ui vanishes but the effective potential sits at the energy of the lowest bound 
state = — I in atomic units. 

Also shown on this plot (black triangles) are numerical results of a Monte-Carlo simu- 
lation of the density matrix. Paths in imaginary time were sampled from the free-particle 
distribution by way of the Levy construction, and the actions of these paths were calculated 
and averaged according to their weights determined from the sampled, free-particle distri- 
bution of paths. See Ceperley [15] for an extremely useful and pedagogical introduction into 
Monte-Carlo methods for calculating A^-particle (Bose) density matrices. 

In fig |2] we show the result for two identical electrons after the correct symmetrization is 
preformed. Again the dotted lines show the result of Ui and the corresponding solids lines 



are the result after anti-symmetrization of the density matrix 26 
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FIG. 1: Effective Coulomb potential for [3 = 1,10. Dotted lines show the Kelbg potential ui{r) 
and corresponding solid line shows potential including correction from lowest bound state. The 
triangles are the results of PIMC (see text). 
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FIG. 2: Dotted lines show the Kelbg potential for a repulsive electron pair. The solid lines show 
the potential after anti-symmetrization of the density matrix. The triangles are the results of the 
PIMC. 



III. EFFECTIVE LATTICE POTENTIALS 



A. Lattice Parameterization 

Now we focus on the problem at hand, obtaining an effective potential appropriate for 
simulations of heavy flavor quarks or plasma quasiparticles. Internal and free energies for 
static quark-antiquark pairs have been calculated at temperatures above deconfinement |8l 
[9], and either of these energies could be used to construct an effective interaction. However, 
it is not clear which potential would best be used in the simulation of dynamical heavy 
quarks. Shuryak and Zahed [ISj have argued that the rapid rotation of the quarks in bound 
states causes the timescales of interest to be relatively short, and therefore the appropriate 
energy for the interaction would be the internal energy. However, Mocsy and Petreczky [T7] . 
have found that potential models which use the free energy instead of the internal energy 
as the interaction potential best describe the spectral functions which were extracted from 
the Euclidean correlators on the lattice using the maximal entropy method. While this does 
not yet rule out the argument of Shuryak and Zahed (there is still significant uncertainty 
in the spectral function obtained using the maximal entropy method), it clearly establishes 
that at the time of this publication, there is no clear answer concerning which potential best 
describes the interaction between heavy quarks. 

With no specific bias towards the use of either energy for describing the interaction, we 
will consider the internal energy. We parameterize the internal energy in the temperature 
range 1.1 — 2Tc as follows: 

(r, T) = _ ^e-'^(^)^ + ^ tanh(a/i(T)r) (28) 
r b 

This potential is similar in form to the screened Cornell potential as first proposed in 
|18j . The first part of the potential is the usual screened Coulomb potential with temper- 
ature dependent coupling (a) and Debye mass (/i). These values were extracted from the 
lattice data and parameterized as follows. For the coupling constant we use a form similar 
to that from the leading log order renormalization group equation, a{T) = [2.4 ln(2T/Tc)]^^. 
The fit to the Debye mass goes as /i oc T which is consistent with the lowest-order 
perturbative calculations [19] with additional corrections to account for the data better: 
/i(r) = 0.0675(T/Te) - (^^)^ 
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FIG. 3: Left: Lattice results for unquenched color singlet internal energy compared to the param- 



eterization in eq. 28 Right: Corresponding QQ potentials from eq. 
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The second part of the potential takes on the same form as a screened hnear confining 
potential. We choose a different functional form than the usual are~^^ in order that the 
potential remains fiat at intermediate and large distances instead of having a maximum and 
then becoming repulsive at larger r as in the case of the Cornell potential. We find that a 
parameterization of the string tension as (T(T)/(GeV^) = 0.32(Tc/T) + 0.8(Tc/T)^® + 0.005 
reproduces the lattice data reasonably well across a large range of temperatures. The final 
parameters are a = 5.5, b = 0.37 GeV and = 0.27 GeV. 

The potential we use is then given as: 



V^i (T, r)= El (T,r)- El (r,oo) 



(29) 



Subtracting away the internal energy at r = oo is convenient when we use the BPL 
convention for the bound state projection operator P' . The results of the fit to the internal 
energy compared to the lattice data as well as the resulting potentials is shown in fig [3| 



B. Effective potentials from temperature-dependent Hamiltonians 

Any parametrization of the lattice results for the heavy quark potentials will be strongly 
dependent on temperature. So far, we have only considered creating effective potentials 
from temperature-independent potentials (for example, the Coulomb potential). Starting 
from temperature-dependent potentials is straightforward. However, because there seems to 



11 



be little discussion of this in the literature, we show how our results may be generalized to 
the temperature-dependent case. 

Consider the following partial differential equation, 

H{f,(3)p{f,7^,T) = -^. (30) 



This equation is identical to the Bloch equation (eq. 12) but with the Hamiltonian, H, 



taken at the fixed temperature /3. Because the potential in if is at a fixed temperature, the 



solution to [30] is the same as before, 

p{f,r,T) = (f\exp{—TH{P)\r). (31) 



By setting t = f3 we find that the right-hand side is exactly the density matrix. 

Therefore the prescription for a temperature-dependent Hamiltonian is simple: fix the 
potential at the temperature of interest for all r, then solve for the density matrix exactly 
as we have done in the previous sections {i.e. either using the Bloch equation or numerically 
with PIMC methods). Finally set r = /3 = 1/T. This example shows how the analogy 
between the matrix elements of a density matrix, and imaginary- time propagators, breaks 
down when the potential has temperature dependence. 

C. Results 

Using the lattice potential from the previous section the methods of Section [TT] can be 
applied. Do the previous results justify using the Kelbg potential to approximate the effective 



QQ potential? Well, since the divergent term in eq. [28j is the screened Coulomb term at 
small r, we may work in units where \Ei{r,T)\ < K The question then is whether or not f3 
is sufficiently small in units of the Bohr energy Eb- The answer is yes, since in the extreme 
case where the Bohr energy is the highest and the temperature is at the lowest (near Tc), 
we actually have T ~ Eb, where our results for the Coulomb case clearly show the Kelbg 
potential working well down to T = O.IEb- So we are confident that our methods will work 
well for the temperature range of interest. 

In fig. [5] we show the results for Ui{r) for charmonium and bottomonium (dotted blue and 
green lines respectively) at a temperature T = 1.2Tc. As in the Coulomb case the effective 
potential takes on a finite value at r = 0. 
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FIG. 4: Lowest QQ bound state energy as a function of temperature. 
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FIG. 5: Effective potential from the lattice data for charmonium and bottomonium at T = 1.2Tc- 
The dotted (upper) curves show ui and the solid (lower) curves show the result after the orthogo- 
nalization to the lowest bound state. 

It is also possible to include the effect of the lowest bound state as was done in the 
Coulomb case. The resultant energy for the lowest bound state of charmonium, bottomo- 
nium and also a light quark bound state having reduced mass M = 0.75 GeV is shown in 
fig. |4j Also shown in fig. [5] as solid lines is the effective potential when the lowest bound 
state is included using the BPL formulation. At T = 1.2Tc the effect of the bound state is 
to lower the potential by about 5% for charmonium and by about 15% for bottomonium. 
The effect is even smaller at higher temperatures. 
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In order for the results given above to be useful in future simulations we now quote 
some useful results for Kelbg-like screened Coulomb and linear potentials. We note that 



equation 19 can be integrated exactly in the case of a screened Coulomb potential, V{r) 



with the result: 



«i(r,A) = ^y^e--e(-^)^/^x [l-e^--2erf(^)+erf(^-^) + e^-erf(^ + ^)] (32) 
and for a linear potential V{r) = br: 



Mr,X) = '1^ X [1 + 3^ - (1 + ^)e-'=/^' + ^(3 + 2^)(l - erf(^))] (33) 



which can be used as an approximation to the confinement term in the lattice potential 
at small distances. However, the inclusion of bound states which are most prominent at low 



temperature and higher masses must be calculated numerically from eq. 22 



IV. CONCLUSIONS 

We have used methods developed for treating effective potentials in electron-ion plasmas 
in order to generate potentials for quarkonium above deconfinement that properly take into 
account quantum effects and symmetry considerations. The nature of the potential changes 
at small distances (r < 1 fm) and reaches a finite value at zero as is also seen in the Coulomb 
case. These effective potentials should be used in classical transport simulations as well as 
in molecular dynamic simulations of the quark-gluon plasma when F is sufficiently small. 
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